Poisson Blending and Mixed-Gradient Blending

Load Image

im_background = imresize(im2double(imread('./data/hiking.jpg')), 0.5, 'bilinear');
im_object1 = imresize(im2double(imread('./data/penguin-chick.jpeg')), 0.5, 'bilinear');
im_object2 = imresize(im2double(imread('./data/penguin.jpg')), 0.5, 'bilinear');

Get a mask on the desired part of the source image

obj1_mask = getMask(im_object1);
Draw polygon around source object in clockwise order, q to stop

Specify a location in the target image where it should be blended.

Transform the source image so that indices of pixels in the source and target regions match.

[im_s1, mask_s1] = alignSource(im_object1, obj1_mask, im_background,0.7);
choose target bottom-center location
im_blend1 = poissonBlend(im_s1, mask_s1, im_background);
figure(3), hold off, imshow(im_blend1);
The above is the result. The penguin blended much more naturally.
obj2_mask = getMask(im_object2);
Draw polygon around source object in clockwise order, q to stop
[im_s2, mask_s2] = alignSource(im_object2, obj2_mask, im_blend1,1);
choose target bottom-center location
im_blend2 = poissonBlend(im_s2, mask_s2, im_blend1);
figure(3), hold off, imshow(im_blend2);
imwrite(im_blend1,'Poisson_blend1.png');
imwrite(im_blend2,'Poisson_blend2.png');
The blending is done quite well. But there are some limitations.

Mixed-Gradient

im_blend_mix1 = mixedBlend(im_s1, mask_s1, im_background);
figure(3), hold off, imshow(im_blend_mix1);
im_blend_mix2 = mixedBlend(im_s2, mask_s2, im_blend_mix1);
figure(3), hold off, imshow(im_blend_mix2);
imwrite(im_blend_mix1,'Mix_blend1.png');
imwrite(im_blend_mix2,'Mix_blend2.png');
function [output] = poissonBlend(im_s,mask_s,im_background)
[imh, imw, nn] = size(im_background);
im2var = zeros(imh, imw);
im2var(1:imh*imw) = 1:imh*imw;
% mask_1의 x,y좌표값 가져오기 -> A의 ncol 계산하기
[y,x] = find(mask_s ~= 0);
[y_s,x_s] = find(mask_s == 0);
% matrix A의 윗부분 만들기
nrow = imh*imw;
ncol = size(x,1)+size(x_s,1);
A = sparse([], [], [], ncol, nrow, 0);
b = sparse(ncol,1);
output = zeros(imh,imw,nn);
for c = 1:nn
t = im_background(:,:,c);
s = im_s(:,:,c);
e = 0;
 
% 마스크 안은 source의 gradient를 따르도록 함.
% 단, 만약 boundary에 걸쳐있다면, target_j의 값과 source_i의 값 비교
for i = 1:size(x,1)
e = e+1;
right_term = 0;
A(e, im2var(y(i),x(i))) = 4;
if mask_s(y(i),x(i)-1) == 0
A(e, im2var(y(i),x(i)-1)) = 0;
right_term = right_term + s(y(i),x(i)) - s(y(i),x(i)-1) + t(y(i),x(i)-1);
else
A(e, im2var(y(i),x(i)-1)) = -1;
right_term = right_term + s(y(i),x(i)) - s(y(i),x(i)-1);
end
 
if mask_s(y(i),x(i)+1) == 0
A(e, im2var(y(i),x(i)+1)) = 0;
right_term = right_term + s(y(i),x(i)) - s(y(i),x(i)+1) + t(y(i),x(i)+1);
else
A(e, im2var(y(i),x(i)+1)) = -1;
right_term = right_term + s(y(i),x(i)) - s(y(i),x(i)+1);
end
if mask_s(y(i)-1,x(i)) == 0
A(e, im2var(y(i)-1,x(i))) = 0;
right_term = right_term + s(y(i),x(i)) - s(y(i)-1,x(i)) + t(y(i)-1,x(i));
else
A(e, im2var(y(i)-1,x(i))) = -1;
right_term = right_term + s(y(i),x(i)) - s(y(i)-1,x(i));
end
 
if mask_s(y(i)+1,x(i)) == 0
A(e, im2var(y(i)+1,x(i))) = 0;
right_term = right_term + s(y(i),x(i)) - s(y(i)+1,x(i)) + t(y(i)+1,x(i));
else
A(e, im2var(y(i)+1,x(i))) = -1;
right_term = right_term + s(y(i),x(i)) - s(y(i)+1,x(i));
end
 
b(e) = right_term;
end
%마스크 밖은 target의 pixel을 따르도록 함
for i = 1:size(x_s,1)
e = e+1;
A(e,im2var(y_s(i),x_s(i))) = 1;
b(e) = t(y_s(i),x_s(i));
end
v = A\b;
n=0;
for j = 1:imw
for i = 1:imh
n = n+1;
output(i,j,c) = v(n);
end
end
end
end
 
function [output] = mixedBlend(im_s, mask_s, im_bg)
[imh, imw, nn] = size(im_bg);
im2var = zeros(imh, imw);
im2var(1:imh*imw) = 1:imh*imw;
% mask_1의 x,y좌표값 가져오기 -> A의 ncol 계산하기
[y,x] = find(mask_s ~= 0);
[y_s,x_s] = find(mask_s == 0);
% matrix A의 윗부분 만들기
nrow = imh*imw;
ncol = size(x,1) + size(x_s,1);
A = sparse([], [], [], ncol, nrow, 0);
b = sparse(ncol,1);
output = zeros(imh,imw,nn);
for c = 1:nn
t = im_bg(:,:,c);
s = im_s(:,:,c);
e = 0;
for i = 1:size(x,1)
e = e+1;
 
right_term = 0;
A(e, im2var(y(i),x(i))) = 4;
 
if mask_s(y(i),x(i)-1) == 0
A(e, im2var(y(i),x(i)-1)) = 0;
right_term = right_term + mix_compare(y(i),x(i),t,s,0,-1) + t(y(i),x(i)-1);
else
A(e, im2var(y(i),x(i)-1)) = -1;
right_term = right_term + mix_compare(y(i),x(i),t,s,0,-1);
end
 
if mask_s(y(i),x(i)+1) == 0
A(e, im2var(y(i),x(i)+1)) = 0;
right_term = right_term + mix_compare(y(i),x(i),t,s,0,1) + t(y(i),x(i)+1);
else
A(e, im2var(y(i),x(i)+1)) = -1;
right_term = right_term + mix_compare(y(i),x(i),t,s,0,1);
end
if mask_s(y(i)-1,x(i)) == 0
A(e, im2var(y(i)-1,x(i))) = 0;
right_term = right_term + mix_compare(y(i),x(i),t,s,-1,0) + t(y(i)-1,x(i));
else
A(e, im2var(y(i)-1,x(i))) = -1;
right_term = right_term + mix_compare(y(i),x(i),t,s,-1,0);
end
 
if mask_s(y(i)+1,x(i)) == 0
A(e, im2var(y(i)+1,x(i))) = 0;
right_term = right_term + mix_compare(y(i),x(i),t,s,1,0) + t(y(i)+1,x(i));
else
A(e, im2var(y(i)+1,x(i))) = -1;
right_term = right_term + mix_compare(y(i),x(i),t,s,1,0);
end
 
b(e) = right_term;
end
for i = 1:size(x_s,1)
e = e+1;
A(e,im2var(y_s(i),x_s(i))) = 1;
b(e) = t(y_s(i),x_s(i));
end
v = A\b;
n=0;
for j = 1:imw
for i = 1:imh
n = n+1;
output(i,j,c) = v(n);
end
end
end
end
 
function [d] = mix_compare(y,x,t,s,a,b)
 
if abs(t(y,x)-t(y+a,x+b)) < abs(s(y,x)-s(y+a,x+b))
d = s(y,x)-s(y+a,x+b);
else
d = t(y,x)-t(y+a,x+b);
end
end
 
function [mask, poly] = getMask(im)
% [mask, poly] = getMask(im)
% Asks user to draw polygon around input image. Provides binary mask of
% polygon and a chain of all interior boundary points.
 
disp('Draw polygon around source object in clockwise order, q to stop')
figure(1), hold off, imagesc(im), axis image;
sx = [];
sy = [];
while 1
figure(1)
[x, y, b] = ginput(1);
if b=='q'
break;
end
sx(end+1) = x;
sy(end+1) = y;
hold on, plot(sx, sy, '*-');
end
 
mask = poly2mask(sx, sy, size(im, 1), size(im, 2));
if nargout>1
[poly.x, poly.y] = mask2chain(mask);
end
end
 
function [im_s2, mask2] = alignSource(im_s, mask, im_t,a)
% im_s2 = alignSource(im_s, mask, im_t)
% Asks user for bottom-center position and outputs an aligned source image.
 
figure(1), hold off, imagesc(im_s), axis image
figure(2), hold off, imagesc(im_t), axis image
[y, x] = find(mask);
y1 = min(y)-1; y2 = max(y)+1; x1 = min(x)-1; x2 = max(x)+1;
im_s2 = zeros(size(im_t));
disp('choose target bottom-center location')
[tx, ty] = ginput(1);
 
yind = (y1:y2);
yind2 = (yind -max(y))*a + round(ty);
yind2 = round(yind2);
xind = (x1:x2);
xind2 = (xind - round(mean(x)))*a + round(tx);
xind2 = round(xind2);
 
y = (y - max(y))*a + round(ty);
x = (x - round(mean(x)))*a + round(tx);
x = round(x);
y = round(y);
ind = y + (x-1)*size(im_t, 1) ;
mask2 = false(size(im_t, 1), size(im_t, 2));
mask2(ind) = true;
im_s2(yind2, xind2, :) = im_s(yind, xind, :);
im_t(repmat(mask2, [1 1 3])) = im_s2(repmat(mask2, [1 1 3]));
 
figure(1), hold off, imagesc(im_s2), axis image;
figure(2), hold off, imagesc(im_t), axis image;
drawnow;
end